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ABSTRACT 


This paper presents the considerations and the resultant approach used to implement design 
sensitivity capability for grids into a large scale, general purpose finite element system (MSC/NASTRAN). 
The design variables are grid perturbations with a rather general linking capability. Moreover, shape and 
sizing variables may be linked together. The design is general enough to facilitate geometric modeling 
techniques for generating design variable linking schemes in an easy and straightforward manner. 

Test cases have been run and validated by comparison with the overall finite difference method. 
The linking of a design sensitivity capability for shape variables in MSC/NASTRAN with an optimizer 
would give designers a powerful, automated tool to carry out practical optimization design of real life, 
complicated structures. 

INTRODUCTION 

This paper presents the considerations and the resultant approach used to implement design 
sensitivity capability for grids into MSC/NASTRAN. MSC/NASTRAN is a large-scale, general purpose 
computer program which solves a wide variety of engineering problems by the finite element method. In 
1983, the design sensitivity analysis (DSA) capability was installed in MSC/NASTRAN. This capability 
has recently been enhanced to include a fully integrated optimization capability for sizing variables. With 
the addition of shape sensitivity capability and with the increasing interest in aerospace and automotive 
industries, this general capability can be used in its own right for carrying out sensitivity analysis of 
complicated real life structures. 

Shape optimization is still in a state where fundamental research is needed (Reference 1). The 
integration of shape optimization concepts within Finite Element Methods (FEM) and Computer Aided 
Designs (CAD) should help to bridge the gap between these two technologies. To be successful, the 
proposed integration of software should lead to a system easy to use. From a practical point of view, the 
computational tool should indeed be employed by design engineers with only a superficial knowledge of 
the theoretical basis of each technique. 
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The following sections will be devoted to the basic procedure for an efficient design variable linking 
scheme using reduced basis concepts. These are necessary to avoid generating unrealistic designs due 
to independent node movements. 

The next section provides an overview of the design sensitivity capability in MSC/NASTRAN. The 
design constraints can be weight, volume, frequency, buckling loads, displacements, stresses, strains or 
forces. The design variables can be shape or sizing variables. The method chosen is a semianalytical 
approach based on a first variation (finite difference scheme) of the system’s equilibrium equations with 
respect to the design variables. 

A brief description follows of the program architecture and considerations that go into implementing 
such a capability into a large-scale, general-purpose computer program. The important considerations 
are ease of use, generality, compatibility with the existing architecture, data organization and manage- 
ment and, finally, a good restart capability for an active man/machine interaction. 

Errors associated with shape sensitivity analysis using the semianalytical approach as shown in 
Reference 2 are examined. An iterative scheme and error index are employed to minimize the errors in 
the sensitivity calculations. 

Two example problems were chosen to validate the capability and to highlight some of the salient 
features. The first example problem is a beam modeled by solid elements, with x-section as design 
variable. The second example problem is a cantilever beam with length as design variable. 

BASIC PROCEDURE 

When dealing with shape optimization problems, the design variables must be selected very 
carefully. The coordinates of the boundary nodes of the finite element model is a straightforward choice. 
This choice, however, exhibits many severe drawbacks. The set of design variables is very large and the 
cost and difficulty of the minimization process increase. It has a tendency to generate unrealistic designs 
due to the independent node movement and additional constraints avoiding such designs are difficult to 
cope with. Moreover, an automatic mesh generator is necessary to maintain the mesh integrity through- 
out the optimization process. One obvious remedy is to avoid a one-to-one correspondence between the 
finite element model and the design variables. 

One way to achieve this goal is to use the concept of "design model" utilizing "reduced basis 
vectors". The general form of this relationship is (Reference 3) 

{Ag} = m {x} (1) 

nxl nxm mxl 

The given design Ag define a system of n variables. We refer to columns of T as basis vectors. 
Clearly the method is most useful if m « n and the method will produce a true optimum only if some 
combination of basis vectors can define that optimum. 
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From Equation 1 , it is obvious that the reduced basis approach is simply a design variable linking 
scheme where the basis vectors are columns of the coefficient matrix T. Ideally, the basis vectors could 
most conveniently be generated using a graphics preprocessor. 

To clarify ideas, consider the coordinate update equation given by 


r^new^ 
| Ynew 

► 

'XouO 

yow 

> + 2 X * ^ 

'D-cosine-x' 

D-cosine-y 

► 

A 

o o 

4—^ 

VZnew^ 


'-Zold> 


^D-cosine-z> 


IcJ 


The constant term vector C x , C y , C z are computed to make necessary adjustments so that x new , 
y new , z new are equal to x old , y oW , z M for the initial values of the design variables x supplied by the user. 
Notice that in Equation 2, for the initial configuration, the values of the initial design variables x supplied 
by the user need not be zero. Equation 2 may be rewritten as 
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or 

{Ag}j = m* {x}j 


this is the same form as shown in Equation 1 . The columns of the T matrix are then the basis 
vectors. The capability exists for the user to input the elements of the T matrix directly or in the future he 
may generate elements of the T matrix using a preprocessor. It should be clear from the above formula- 
tion that the design variables are the grid perturbations and not the grid coordinates. 

DESIGN SENSITIVITY CAPABILITY IN MSC/NASTRAN 

Design sensitivity analysis estimates the effects of interrelated design variables such as element 
properties and materials on the structural response quantities such as displacement, stress, natural 
frequency, buckling loads - and for composites lamina stresses and failure indices. Design sensitivity 
coefficients are defined as the gradients of the design constraints with respect to the design variables at 
the current design point. The method chosen for incorporation into MSC/NASTRAN is a semianalytical 
approach, based on a first variation (finite difference scheme) of the systems equilibrium equations with 
respect to the design variables. 

Let Vi (bj , u) be a set of design constraints which are functions of bj design variables and dis- 
placements u. The design constraints are expressed as 


Vi (bj , u) < 0 


The first variation in v is given as 



consider u as a function of bj, then 


(8u} - [!;]<*» 
nxj jxl 


and, therefore, 


or 
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The matrix 55- can be evaluated by taking the first variation of the systems equilibrium equation, 

dbj 


[K] {u } = {P } 


which gives 


[K ] {Au } +[AK]{u} ={AP} 


solving for {Au } 


{Au } = [K J' 1 ({AP } - [AK ] {u }) 


or 
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[Au ] = [K r 1 ({AP(Ab 1 ) } , (AP(Ab2) } {AP(Abj) }) 


- [K r 1 ([AK(Ab,) ] {U } , [AK(Ab 2 ) ] {U } [AK(Abj) ] {U }) 

The elements of matrix for an element constraint such as stress, force, or failure index can 

be expressed by the relationship 


{Vi} »[S] e {u} 


or 

[$]-»» 


The design sensitivity coefficient matrices may thus be expressed as 


M ■({£}•{£> {$})i- 

♦m({£M£} 


From this equation it is easy to see that the number of additional case control records (additional 
loading cases) required for design sensitivity analysis is equal to the number of design variables for each 
subcase (Design Space Technique). 

A typical term of the coefficient matrix may be written as 


( s b+ab u b 
AB 


S B u B \ /S B u B+AB 
AB / + V AB 


S B u B \ 
AB / 


where B represents the base line or original state and B + AB represents the perturbed state. The first 
expression in parentheses on the right-hand side is thus the change in response quantity due to a 
change in design variable for the original solution vector. The second term represents the change in 
response quantity due to a change in displacement for the unperturbed design variable. For displace- 
ment constraints, the first term in parentheses on the right-hand side is identically zero. 

The design constraints can be weight, volume, frequency, buckling loads, displacements, stresses, 
strains, forces, lamina stresses, lamina strains, failure indices or user supplied synthetic equations. The 
design variables can be grid movements or properties. The shape and sizing variables may be linked 
together. 


702 



PROGRAM ARCHITECTURE 


In order to understand the reasons behind how a development is introduced into a large finite 
element program, a knowledge of the program architecture and technical purpose is necessary. A brief 
description of MSC/NASTRAN is presented as background (Reference 4) 

The cornerstone of MSC/NASTRAN’s architecture is its Executive System, whose essential 
functions are to establish and control the sequence of calculations, to allocate files, and to maintain a 
restart capability. Engineering calculations are performed by approximately 200 Functional Modules 
which communicate only with the Executive System and not with each other. Flexibility is maintained by 
a macro-instruction language called DMAP, which is under user control, but which also serves to estab- 
lish preformatted calculation sequences for the major types of analysis, including linear analysis, buck- 
ling, vibration mode analysis, and design sensitivity. 

The calculation of finite element data is concentrated exclusively in a few modules. The element 
matrices for stiffness, structural damping, and differential stiffness for elements of the structural model 
are generated in the Element Matrix Generator (EMG) module. These element matrices are subse- 
quently assembled to form the elastic stiffness matrix, the structural damping matrix, the mass matrix, or 
the differential stiffness matrix. 

The element contribution to the load vector is generated in the load generator module and the 
element stress and force are generated in the recovery module. In all these modules, the finite element 
descriptions are defined in the Element Summary Table (EST). The EST contains the element connec- 
tion, material property and sectional property information. 

Taking advantage of the table driven concept used by the element modules, much of the element 
dependent development can be avoided in implementing design sensitivity. The reason is that a proce- 
dure could be developed which only involves building EST tables that would cause existing modules to 
form the necessary element data. 

How a given capability is introduced into a commercial general purpose finite element program is 
as important an issue to the user as its theoretical sophistication. If the user views a capability as hard to 
use, as having an insufficient capacity to solve his problem, or taking an inordinate time to comprehend 
its output, the product is of little practical use. In addition, the program developer, while heeding the 
user’s needs, has to keep sight of the program as a whole when adding new capabilities. This involves 
interfacing well with existing capabilities, maintaining program reliability and generality, and producing 
software that makes effective use of computer resources. 
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The user interface is a major consideration in the design of a new capability. The following issues 
were considered when building up the design sensitivity analysis (DSA) capability. 

1. DSA input should be straightforward, but allow flexibility to model complex structural design 
concepts. 

2. DSA output should be concise and easily understood. 

3. Avoid arbitrary program limits which restrict the allowable element types, constraint quantities, 
and problem size. 

4. Provide an interface for external optimization postprocessors. 

A brief discussion of the processes involved in a typical DSA STATIC analysis in MSC/NASTRAN 
will help bring into perspective the work involved in the various parts of the DSA solution. 

DSA in a STATICS analysis is based on solving for {Au} the first order variation of the nodal 
equilibrium equation: 


[K°] {Au } = {AP } - [AK] {u 0 } 

The DSA problem in this paper is considered to be the additional task required after the solution of 
primary analysis. By restarting from the primary STATIC analysis, the solution of the DSA system 
equation only involves the calculation of the right-hand side and the backward pass operation in the 
solution of Au. 

The woik involved in solving the system equations (backward pass operation) is a function of the 
product of the number of design variables and loading conditions. The following DSA tasks are required 
in addition to solving the system equations: 

1 . DSA Data Organization 

2. DSA Data Assembly 

3. DSA Data Recovery 

These tasks are functions of the triple product of the number of design variables, design constraints and 
loading conditions. For large DSA problems, the data organization, assembly and recovery tasks are the 
dominant users of computer resources. 

Another major consideration was to support all structural finite element types in MSC/NASTRAN. 
Since a large number of the elements developed are semiempirical, the determination of consistent 
element derivative formulations cannot be practically accomplished. Therefore, a method was developed 
to calculate element derivatives by a differencing scheme about the current design point. This method 
involved the calculation of the element matrix at the design point plus or minus the user specified design 
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variable increment. This element data is differenced with the data at the design point to determine the 
corresponding element derivatives. For example, the following shows the change in element stiffness 
due to a change in the design variable. 




An initial analysis is carried out to identify critical constraints and a data base is created. In the 
succeeding run, information about constraints, design variables, maximum and minimum side constraints 
is supplied. A special OMAP package was created which exploits the data base technology. 

ERRORS ASSOCIATED WITH SEMI-ANALYTICAL APPROACH IN SHAPE SENSITIVITY 

In Reference 2 it was shown that the semi-analytical approach can have serious accuracy prob- 
lems for shape design variables in structures modeled by beam, plate, truss, frame and solid elements. 
An error index was developed to test the accuracy of the semi-analytical approach and some methods 
were proposed for improving the accuracy of the semi-analytical method. In the following section, the 
interactive scheme proposed in Reference 2 is examined in greater detail. 

Consider the systems equilibrium equation 

[k]{u} ={P} 

Taking the first variation of the systems equilibrium equation gives 

[k + Ak] {u + Au} - {P + AP} 

Expanding the above equation and retaining the second order terms gives 

[k] {Au} * {AP} - [Ak] {u} - [Ak] {AU} 
which can be cast as an iterative scheme 

[k] (Au}j = {AP} - [Ak] {u} 0 - [Ak] {Au} m (3) 

This iterative scheme can be used quite effectively with the error index as shown in Reference 3, 
i.e., 


^imn _ 


[A] 1 - [A] m 
[A]" 


(4) 
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An error index of the type described above is almost a must for a large scale system in the op- 
timization context. If the sensitivity derivatives are significantly in error, the program needs to detect it 
and stop the execution to save time and resources of the user. 

NUMERICAL EXAMPLES 

Two example problems were chosen to validate the capability and to highlight some of the salient 
features. 

Example 1 - Beam using solid elements 

A cantilever beam is subjected to a tip loading. The model of the beam is shown in Figure 1 . The 
beam is modeled using solid hexahedron elements. The analysis model consists of five HEXA elements 
and 24 grids. 
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The design model consists of one design variable to perturb beam cross section to maximize the 
bending inertia about the y-axis. The cross-section is perturbed as shown in Figure 2 below. 



Figure 2. Design Variable Linking Scheme. 

P y , and P Z| are perturbation vectors in the y and z direction respectively. 

The perturbation vectors are given by 

P Yl =Py 2 = O.Ol(o? + j + Ok) 

Py 3 =Py 4 =0.0l(0i-1+0k) 

P*i =Pz 3 =O.Ol(oi +ol + k) 

Pz 2 =Pz 4 =O.Ol(oi+Oj-k) 

The sensitivity coefficient results calculated using the semi-analytical approach are compared to 
the Overall Finite Difference (OFD) approach, wherein the entire problem is solved again for the per- 
turbed configuration. The results are shown in Table 1. As can be seen, the correlation between he 
semi-analytical approach and the OFD is excellent for this particular example problem. 
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Table 1. Vertical displacements along beam length.. 

Example 2 • Cantilever beam subjected to end moment 

The same cantilever beam that was considered in Reference 2 is taken as an example here. The 
cantilever beam has uniform rigidity El and length L under a tip moment M as shown in Figure 3. 




1 


» •- 

Id) (2) 
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Figure 3. Cantilever beam subjected to end moment. 


The modulus of elasticity is taken equal to E - 10 7 . The theoretical tip deflection for the above 
configuration is 


5 = 


Mi 2 

2EI 


(5) 


The structure has been idealized into 20 beam elements. Sensitivity coefficients have been 
calculated for the tip displacement with respect to the length of the beam. The design model contains a 
single design variable, i.e., the length of the beam. The grid perturbations of all the grids have been 
linked together so that the perturbations vary linearly from roof to tip of the beam. 
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Parametric studies have been carried out to determine the effect of the step size on the error. 
From Equation 5, the exact results are available. 

First the sensitivity analysis is carried out using the Overall Finite Difference (OFD) and the 
semi-analytical (SA) method for step sizes of 1%, 0.5%, 0.1 % and 0.01%. The results, without using any 
iterative schemes are shown in Table 2. As can be seen, the results are quite accurate for a step size of 
0.01%. Whereas they progressively degrade for increasing values of the step size and become quite 
unacceptable for a step size of 1%. 

In Reference 2. the errors for the beam-type structure are associated with an incompatiblity of the 
sensitivity field with the structural model. The error in the finite difference approximation consists of a 
truncation error because of neglecting higher order terms in the Taylor series expansion and a condition 
error because of the limited precision available for the computer. Thus, an optimum value of step-size 
would minimize the truncation error without the condition error becoming significant. As suggested in 
Reference 3, central difference scheme is an alternative to iterations. 


Step-size (%) 

Sensitivity Coefficient 
A 

1.0 

-10.459 

0.5 

-4.716 

0.1 

0.010 

0.01 

1.08 

0.001 

1.1885 

0.0001 

1.1874 

0.00001 

0.914* 


* Degrades. (Theoretical value = 1 .2) 

Table 2. Variation of sensitivity coefficients with respect to step-size. 


Next, we use the iterative scheme of Equation 3 to converge to the correct solution. As can be 
seen from Figure 4, the smaller the value of the step size, the lesser the number of iterations required to 
converge to the correct solution. 
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Figure 4. Plot of number of iterations versus error in tip-displacement using SA method. 


CONCLUSIONS 

This paper presented the considerations and the resultant approach used to implement design 
sensitivity capability for grids into a large scale, general purpose finite element system (MSC/NASTRAN). 
The design variables are grid perturbations with a rather general linking capability. Moreover shape and 
sizing variables may be linked together. The design is general enough to facilitate geometric modeling 
techniques for generating design variable linking schemes in an easy and straightforward manner. 

The errors shown to be associated with the semianalytic method for shape variables for beam type 
structures can be mitigated by resorting to an iterative scheme. Examples have been presented high- 
lighting the salient features of the approach. 
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